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Animals of the same species exhibit similar behaviours that are advantageously 
adapted to their body and environment. These behaviours are shaped at the species 
level by selection pressures over evolutionary timescales. Yet, it remains unclear 


how these common behavioural adaptations emerge from the idiosyncratic neural 
circuitry of each individual. The overall organization of neural circuits is preserved 
across individuals’ because of their common evolutionarily specified developmental 
programme? *. Such organization at the circuit level may constrain neural activity> $, 


leading to low-dimensional latent dynamics across the neural population 


9-11 


Accordingly, here we suggested that the shared circuit-level constraints within a 
species would lead to suitably preserved latent dynamics across individuals. We 
analysed recordings of neural populations from monkey and mouse motor cortex 

to demonstrate that neural dynamics in individuals from the same species are 
surprisingly preserved when they perform similar behaviour. Neural population 
dynamics were also preserved when animals consciously planned future movements 
without overt behaviour” and enabled the decoding of planned and ongoing movement 
across different individuals. Furthermore, we found that preserved neural dynamics 
extend beyond cortical regions to the dorsal striatum, an evolutionarily older 
structure”, Finally, we used neural network models to demonstrate that behavioural 
similarity is necessary but not sufficient for this preservation. We posit that these 
emergent dynamics result from evolutionary constraints on brain development and 
thus reflect fundamental properties of the neural basis of behaviour. 


The behaviour of each individual in a species is driven by the coordi- 
nated activity of neural populations throughout the brain. This activity 
emerges from the latent dynamics, which are the time-dependent acti- 
vation of the dominant patterns of neural covariation’”. These latent 
dynamics seem to be shaped by circuit and biophysical constraints ®. 
Given the large differences in brain circuitry across individuals from the 
same species—including in some type-specific neurons, dendritic mor- 
phology and receptor distribution °—it remains unclear how similar 
adaptive behaviours emerge from such idiosyncratic neural circuitry. 
One possibility is that unique circuits in each individual generate unique 
latent dynamics that produce the same behavioural output. Indeed, 
the high degrees of freedom of neural activity relative to behaviour?” 
could allow distinct latent dynamics to produce similar behaviour. 
Alternatively, the same behaviour performed by two individuals could 
be produced by preserved latent dynamics. This preservation would 
emerge from the common organization of neural circuits across indi- 
viduals resulting from a species-specific developmental programme. 

Here, we adopt the last hypothesis: different individuals from the 
same species engaged in the same behaviour generate preserved 


latent dynamics. We posit that preserved circuit constraints give 
rise to a species-wide neural landscape and the individual-specific 
latent dynamics observed during a behaviour are different instan- 
tiations ofa common trajectory through this landscape (Fig. 1). Our 
hypothesis provides several testable predictions. First, because 
low-level details of neural circuits are idiosyncratic, they should not 
be necessary to account for the emergence of species-typical behav- 
iours. Accordingly, different animals of the same species engaged 
in the same behaviour should exhibit preserved latent dynamics. 
Second, the extent of preservation of the latent dynamics across 
individuals should be constrained by the similarity of the behav- 
ioural output. Third, because low-dimensional latent dynamics 
are found throughout the brain, not just in cortical regions” *, we 
should also observe preserved latent dynamics in structures that 
have co-evolved with cortex for hundreds of millions of years such 
as the basal ganglia. Fourth, because covert behaviours seem to be 
mediated by the same neural circuits as overt behaviours”, we should 
find shared latent dynamics across animals performing the same 
cognitive task. 
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Fig. 1| Hypothesis. a,b, Different individuals from the same species 
performing the same behaviour will generate preserved neural population 
latent dynamics by instantiating a species-wide ‘neural landscape’ embedded 
in D dimensions (Dim) of neural activity. c, These preserved latent dynamics 
can be revealed by ‘aligning’ the latent dynamics estimated from neural 
population recordings of each individual (N,,-,and Ng, 3 illustrate three 
neurons recorded from monkey A and monkey B, respectively). 


Preserved latent dynamics across animals 


We tested the four predictions outlined above using neural popula- 
tion recordings from monkeys and mice as they performed upper 
limb tasks. First, we analysed motor cortical recordings from mon- 
keys engaged in an instructed-delay centre-out reaching task with 
eight targets (Fig. 2a,b and Extended Data Fig. 1a; Methods). All three 
monkeys were well-trained in the task and exhibited highly stereo- 
typed hand trajectories (mean trajectory correlation for each monkey: 
r=0.89, 0.90 and 0.92; Extended Data Fig. 1b). For each session, we used 
principal component analysis (PCA) to estimate the latent dynamics 
underlying overt movement execution by projecting the firing rates of 
each recorded neuron (or multi-unit, for monkey J) onto the leading ten 
PCA axes (the neural modes; examples in Fig. 2c; Extended Data Fig. 1c). 
We then aligned the latent dynamics of each pair of experimental ses- 
sions from two different animals using canonical correlation analysis 
(CCA?*),a method that maximizes the correlations between two sets of 
signals through linear transformations (similar to refs. 26-28). 

This linear method revealed that the ostensibly different latent 
dynamics of two different monkeys are indeed highly preserved 
(Fig. 2e). The across-animal correlations greatly exceeded two lower- 
bound controls. The first was established by aligning randomly selected 
behavioural epochs (‘control’ in Fig. 2f) and the second was based on 
surrogate data that conserved the statistical structure of the neural 
activity” (tensor maximum entropy (TME); Extended Data Fig. 2a,b). 
Most importantly, these correlations were nearly as high as values 
obtained by aligning two subsets of trials within a single session from 
the same animal (‘within’ in Fig. 2f; further examples in Extended Data 
Fig. 3a). This result further held across all pairs of sessions from all 
three monkeys (n = 126 sessions; Fig. 2g). The aligned neural modes 
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captured a large fraction of neural variance (Extended Data Fig. 2e) and 
the results did not depend onthe assumed dimensionality of the neural 
manifold (Extended Data Fig. 2f) or the alignment method (Extended 
Data Fig. 2g,h). 

Although we have demonstrated the preservation of latent dynam- 
ics across animals, these shared dynamics may not necessarily be 
relevant to behaviour. To address this, we trained neural network 
decoders” (long short-term memory networks (LSTMs)) to predict 
the hand trajectories of one animal and tested their performance ona 
second animal. The across-animal decoding accuracy approached the 
upper bound provided by the performance of decoders trained and 
tested on the same animal (Fig. 2h; example predictions in Extended 
Data Fig. 4). Thus, the preserved motor cortical dynamics across 
animals contain detailed information about the ongoing movement 
kinematics. 

We then analysed data from four mice trained to perform a reaching 
and pulling joystick task (Fig. 2i and Extended Data Fig. 5). We found 
that both within and across individuals, the behavioural output was 
less similar from trial to trial than for the monkey dataset (compare 
Extended Data Figs. 1a,b to 5a,b). As predicted, our alignment proce- 
dure revealed that motor cortical latent dynamics were largely pre- 
served across mice both when reaching to two different targets (Fig. 2j 
and Extended Data Figs. 3b and 2c; example in Extended Data Fig. 3e) 
and when subsequently pulling at two different force levels (Extended 
Data Fig. 6a). Yet, these correlations were lower than those of monkeys 
(compare Fig. 2g to 2j), which directly impacted our ability to decode 
movement kinematics across animals (Extended Data Fig. 6b). This 
difference could be explained by the more highly stereotyped behav- 
iour of the monkeys compared to the mice (Fig. 2k inset). Comparing 
between species confirmed that behavioural stereotypy was associated 
with both the preservation of the latent dynamics across individuals 
(Fig. 2k) and across-animal decoding accuracy (Extended Data Fig. 6f). 
These results demonstrate in two evolutionarily divergent species that 
there is a direct correspondence between the similarity of behavioural 
output and the preservation of motor cortical latent dynamics across 
individuals. 


Necessity of behavioural similarity 


In the preceding analyses, we studied tasks comprising a few condi- 
tions that inadvertently imposed a topological structure in the pro- 
duced movements”. We sought to establish that preserved latent 
dynamics do not merely reflect this structure. First, we tested whether 
this topological structure is sufficient to produce preserved latent 
dynamics and found that the preservation across individuals was 
significantly impaired (Extended Data Fig. 7). We then further dem- 
onstrated the existence of preserved latent dynamics during a con- 
tinuous and less-structured task in which monkeys rapidly generated 
sequences of random reaches?!” (Fig. 3a). The produced movements 
were highly varied with little organization (Fig. 3b and Extended Data 
Fig. 8a). To facilitate alignment, we parcellated the workspace to match 
movements across individuals on the basis of the initial hand position 
and reach direction to generate as many as 29 similar reaching condi- 
tions (Fig. 3c). Despite the dramatic increase in number of conditions 
and behavioural variability, we could still uncover preserved latent 
dynamics across individuals (Fig. 3d,e and Extended Data Fig. 2d). The 
increased complexity of behavioural output in this task allowed us to 
directly study the relationship between the number of conditions and 
the preservation of latent dynamics. We subsampled the conditions 
and found that latent dynamics across individuals were preserved 
for the entire range (up to 29) considered (Fig. 3f and Extended Data 
Fig. 8b) even when as few as about 20 neurons are included (Fig. 3g and 
Extended Data Fig. 8c). The preservation of latent dynamics across 
individuals decreased when we shuffled the conditions (Extended 
Data Fig. 8d), thereby reducing behavioural similarity. This agrees with 


Behavioural task and neural recordings in monkeys 
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Fig. 2| Preserved latent dynamics across animals performing the same 
behaviour. a, Monkeys performed an eight-target instructed-delay centre-out 
reaching task using a planar manipulandum. b, Example normalized neural 
firing rates aligned to movement onset (top) and hand trajectories (bottom) 
for two monkeys. Each column shows one reach to the eight targets indicated 
by the arrows. Note that monkey C received two sets of motor cortical implants, 
C, and Cp, with C, denoting the implant in the left hemisphere and Cgthe implant 
intheright hemisphere. c-e, Three-dimensional representation of the motor 
cortical latent dynamics for the two monkeys plotted in b before (c) and after 
(e) alignment with CCA (d). f, Correlations of the aligned (red) latent dynamics 
for the example comparison in c-e compared to within-monkey correlations 
(grey) and alower-bound control (orange). g, Preserved latent dynamics across 
all pairs of21 sessions (n=126 comparisons) from three different monkeys. 
Histograms show the mean correlation across the leading four dimensions. 
Line and shaded area, mean + s.d. h, Decoders trained on aligned latent 
dynamics from one monkey predict continuous hand kinematics of a different 


monkey (blue). Results compared to decoders trained and tested within the 
same session (grey) and without alignment (green). Data points, individual 
comparisons (n = 126) between sessions from different monkeys. Error bars, 
mean +s.d. Statistical tests: two-sided Wilcoxon’s rank sum test, P=3.1 10° 
for decoding performance between across-animal correlation and within- 
animal correlation, P=2.0 x 10” between across-animal correlation and lower 
bound. i, Mice grasped and pulled a joystick in two positions (left or right). 

j, Preserved latent dynamics across mice performing the grasping and pulling 
task. Data include six sessions across four different mice (n=13 comparisons), 
formatted asing. k, Preserved latent dynamics across animals is related to the 
similarity of their behaviour. For each session pair, the mean of the top four 
canonical correlations (CCs) between the latent dynamics against the mean 
behavioural correlation. Single dots, pairs of sessions colour-coded by species. 
Inset: behavioural correlations for all pairs of trials from different mice and 
monkeys. Circles, mean. 


Nature | www.nature.com | 3 


Article 


Identifying preserved latent dynamics during sequential reaches 
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Preserved latent dynamics during matched reaches to random targets 


d Example session e 
= Across monkeys == Across monkeys 
1.0 Within monkey C, 1.0 Within monkey 
Within monkey M == Control 
c 087%  — Control 5 0.8 
fe} N 2 
ğ 06 5 06 — 
xo) g 
5 04 5 04 
oO O 
0.2 me a aa D 
0 (0) 
a S S. ee | 
12 4 6 8 10 12 4 6 8 10 


Neural mode Neural mode 


Preservation of latent dynamics is minimally dependent on the number of 
conditions or neurons 


f io EEOAE g 1.0 
xample session Example session 
g£ os \ se os . \ 
os os 
25 os i E] rE 
S 04 8S 04 
§ 5 § 
$38 02 $8 02 
0 0 
T T T T T f Į js T 
3 9 15 21 27 10 50 90 130 


Number of conditions Number of neurons 


Fig. 3 | Preserved latent dynamics in a sequential reaching task with higher 
behavioural complexity. a, Monkeys were trained to perform sequences of 
four reaches to randomly placed target locations using the planar manipulandum. 
b, Example hand positions in the workspace for all reaching movements made 
during a session from two example monkeys. c, We matched similar movements 
between the two monkeys on the basis of their starting location, duration and 
direction. Each colour represents a different representative condition and 
matched reaching trajectories are overlaid for the two monkeys. Circles, 
starting position. d, Preserved motor cortical latent dynamics across monkeys 
performing matched reaches for the example sessions inc. Note that the 
correlation between aligned latent dynamics across monkeys was quite close 
to the within-monkey correlations (grey) and largely exceeded our lower- 
bound control (orange). e, Preserved motor cortical latent dynamics across all 
pairs of sessions from three different monkeys. Data include n =10 comparisons 
and are presented as in Fig. 2g. Line and shaded area, mean + s.d. f, Preserved 
latent dynamics hold for a wide range of conditions. The plot shows, for each 
pair of sessions, the mean of the top four CCs between the latent dynamics as 

a function of the number of conditions subsampled from the total available in 
each session. Dark trace, example session from b-d. g, Preserved latent 
dynamics can be uncovered from even sparsely sampled neural populations. 
The plot shows, for each pair of sessions, the mean of the top four CCs between 
the latent dynamics as a function of the number of neurons subsampled from 
the total available in each session. Dark trace, example session from b-d. 


the previous comparison of mice and monkeys: monkeys performing 
amore complex task (eight conditions) than mice (two to four condi- 
tions) had higher preservation of latent dynamics because of the more 
stereotyped behavioural output. 

We performed more control analyses to confirm that this rela- 
tionship between behavioural stereotypy and preservation of latent 
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dynamics across individuals fromthe same species is not a trivial con- 
sequence of our methodology. First, we compared our results with a 
previous study that investigated the dynamics of motor cortical activity 
within an individual across different but related wrist manipulation 
and reach-to-grasp tasks” (Extended Data Fig. 9). During each task, 
individuals activated the same muscles in a slightly different man- 
ner (Extended Data Fig. 9b-d). Accordingly, the latent dynamics of 
the same monkey performing two distinct but related behaviours were 
much less preserved than those of different monkeys performing the 
same behaviour (Extended Data Fig. 9g). Last, as we have previously 
shown, preserved latent dynamics cannot be explained by stable move- 
ment tuning”, nor do they persist following nonlinear transforma- 
tions”®. These analyses demonstrate that the alignment method alone 
is not sufficient to uncover preserved latent dynamics, even within 
individuals. 


Preserved dynamics in dorsal striatum 


Given that the motor cortex is the main cortical output to the spinal 
circuits that generate movement, the close correspondence between 
motor cortical latent dynamics and behavioural output may uniquely 
result from the architecture and projections of this region. To test 
whether preserved latent dynamics exist across the brain, we studied 
the subcortical nuclei of basal ganglia, which do not directly project 
to spinal cord but are crucial for various aspects of behaviour** **. We 
predicted that basal ganglia latent dynamics would also be preserved 
across animals performing the same task. Replicating our alignment 
analysis on neural population recordings from mouse dorsolateral 
striatum during a reaching and pulling task (Fig. 4a) revealed preserved 
latent dynamics across individuals (Fig. 4b and Extended Data Fig. 3c; 
example in Extended Data Fig. 3f). 

Moreover, despite the vast differences in circuit and cellular archi- 
tecture between motor cortex and striatum**”’, both the across- 
animal correlations (compare Fig. 4b and Fig. 2j) and the across- 
animal decoding performance of hand trajectories (Fig. 4c and 
Extended Data Fig. 6f) were equally large for both regions (Extended 
Data Fig. 6c-e). Therefore, stable latent dynamics across ani- 
mals performing the same behaviour are not confined to motor 
cortex—they extend to different regions throughout the brain, 
including an evolutionarily older structure that is shared among all 
vertebrates”. 


Preserved dynamics during covert behaviour 


We have shown the preservation of latent dynamics across brain 
regions during active, overt behaviour. However, animals also engage 
ina variety of covert behaviours such as deliberation and planning. 
These processes require neural activity that is predominantly internally 
generated by the brain. Given that such covert behaviours involve brain 
regions that also mediate overt behaviours”>”°, we predicted that the 
latent dynamics underlying these more cognitive processes would 
also be shared across individuals of the same species. We tested this 
prediction by analysing motor cortical activity as monkeys planned an 
upcoming movement before executing it (Fig. 4d). The latent dynam- 
ics during the instructed-delay period were highly correlated across 
animals (Fig. 4e and Extended Data Fig. 3d; example in Extended Data 
Fig. 3g) and were virtually identical to those during overt reaching 
behaviour (compare Fig. 2g). Moreover, these aligned latent dynam- 
ics were also predictive of behaviour: Bayesian models predicting 
the upcoming reaching target based on the aligned latent dynamics 
from one monkey generalized to another (Fig. 4f). Thus, different 
individuals use preserved latent dynamics not only to execute the 
same movement but also to perform the same covert mental process. 
This result also strengthens our previous observation of preserved 
latent dynamics during overt behaviour. Afferent feedback arriving 
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Preserved latent dynamics during 
covert motor planning 
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Fig. 4| Preserved latent dynamics in the basal ganglia and during covert 
behaviour. a, We studied whether shared latent dynamics can also be found 

in subcortical regions by analysing recordings from mouse dorsolateral 
striatum during a reaching and pulling task. b, Preserved striatal latent 
dynamics across mice performing the same grasping and pulling task. Data 
include n =13 comparisons across four mice and are presented as in Fig. 2g. Line 
and shaded area, mean + s.d. c, Hand trajectories can be decoded across mice 
using their preserved striatal latent dynamics. Data presented as in Fig. 2g. 
Error bars, mean + s.d. Statistical tests: two-sided Wilcoxon’s rank sum test, 
P=3.0 x 10-8 for decoding performance between across-animal correlation 

and within-animal correlation, P= 7.3 x 10™* between across-animal correlation 
and lower bound. d, To investigate the preservation of latent dynamics during 
covert behaviour, we examined the preparatory period preceding movement 
execution in the monkey centre-out task. e, Preserved motor cortical latent 
dynamics when preparing to reach toa target. Data include 18 sessions and 
n=72comparisons across two different monkeys and are presented asin 

Fig. 2g. Line and shaded area, mean + s.d. f, Preserved latent dynamics during 
covert behaviour contain behaviourally relevant information. Naive Bayes 
classifiers trained on the aligned latent dynamics (blue) predicted the intended 
target ina different monkey virtually as well as classifiers trained and tested on 
the same monkey (grey), whereas classification was poor without alignment 
(green). Error bars, mean + s.d. Statistical tests: two-sided Wilcoxon’s rank sum 
test, P=1.6 x 10‘ for classifier performance between across-animal correlation 
and within-animal correlation, P=1.7 x 10” between across-animal correlation 
and lower bound. 


at the motor cortex*!” could partially explain the observed similarity 
in latent dynamics during overt movement, yet the latent dynamics 
are entirely internally generated during covert processes such as 
movement planning. 


Behavioural similarity is not sufficient 

Our hypothesis requires that behavioural similarity is necessary but 
not sufficient to allow for alignment of latent dynamics. To test this, 
we trained recurrent neural networks (RNNs) to produce similar 


behavioural output while generating distinct latent dynamics. We 
devised an RNN simulation in which we had control over the degree 
of preservation of latent dynamics by varying a parameter of the cost 
function, a (Fig. 5a and Extended Data Fig. 10). We thus created pairs 
of models that generated highly similar behaviour (Fig. Sb-d) but 
exhibited distinct latent dynamics as evidenced by the relative lack of 
alignment (Fig. 5e). We predicted that this decrease in preservation of 
latent dynamics would be driven by differences in underlying circuit 
properties. When we reverse-engineered the weights of the different 
networks, we found that more dissimilar latent dynamics corresponded 
to larger changes in the variance and dimensionality of the weights 
changes during training (Fig. 5f). Thus, preservation of latent dynam- 
ics is not just a trivial consequence of behavioural similarity; instead, 
it probably reflects fundamental organization and constraints in the 
underlying circuit implementation. 


Discussion 


Neural population latent dynamics have been proposed as first-level 
explainers of behavioural and cognitive functions””, a view that has 
shed light onto the neural basis of numerous phenomena, suchas pro- 
cesses underlying covert? and overt behaviour®**“5, how learn- 
ing may happen in neural circuits**°*” and how information may flow 
between different brain regions*”**. Here, we extend recent works??? 
to show that latent dynamics are shared across different individuals 
engaged in the same behaviour for a range of behavioural complexity. 
The discovery of preserved latent dynamics across individuals will 
impact both fundamental and applied neuroscience, in particular the 
development of brain-controlled devices such as neuroprosthetics*™: 
with proper alignment, decoders trained on one participant could be 
readily translated to other individuals *’ to minimize training and 
deployment time. 

We studied two regions (motor cortex and dorsolateral striatum) 
whose functions are tightly linked to the production of limb move- 
ments. In the behaviours studied here, we demonstrated preserved 
latent dynamics using a relatively few dimensions. This low dimension- 
ality could be attributed in part to the temporally smooth and relatively 
constrained movements in these datasets. Further, the limbs of dif- 
ferent individuals from the same species share similar biomechanical 
properties that, throughout evolution and development, have imposed 
extra constraints enforcing the preservation of latent dynamics. By 
contrast, highly sensory-driven regions with rapidly changing inputs 
may require a higher level of granularity in the analysis but we none- 
theless expect to see similar preservation of latent dynamics across 
individuals. ‘Higher’ brain regions (for example, frontal cortex) that 
serve more abstract cognitive functions may show less preservation of 
latent dynamics across animals because of differences in their internal 
states (attention, motivation, satiation and so on). Even if the most 
task-relevant aspects of the latent dynamics were preserved across 
individuals, the influence of these more abstract internal states“ 
could alter the overall latent dynamics, making the activity of these 
regions less amenable to alignment. Moreover, different animals may 
also use different covert strategies to solve the same cognitive task 
based on their biases acquired through learning and past experience™. 
Yet, we expect that the neural population latent dynamics would still 
be preserved if it were possible to appropriately match the internal 
states or strategies across individuals. 

The influence of learning and experience on neural circuit organiza- 
tion and the resulting latent dynamics is a fascinating open question. 
Throughout an individual’s lifetime, new skills are acquired through 
practice, which are consolidated through changes in neural circuitry. 
what, if any, differences should exist in the latent dynamics produced 
by the motor cortices of, say, a virtuoso guitarist and a first-year guitar 
student? Regardless of whether those circuits were shaped by practice 
or even development, we propose that neural circuits are tuned to 
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Training ‘constrained’ neural networks to establish that similar behaviour 
is not sufficient for the preservation of latent dynamics 


a b 
Step 1: train standard network 7 lee a 
‘Standard 
Target cue _, : N £ i Ap ; 
eocccecce Ces > a 
x, y) = 
Hold cue > eeu 
E 
fi | 
Step 2: train constrained network w 
f] 
o AWI 
Target cue œ ) o 
eoeeeees es ot 
Hold cue > Ə f ay a = 0.50 


Constrained networks generate similar behaviour output 


s 0.3 d 1.0 
= E 
a 5.9 
5; o% È ee r à 
= 0.2 >p 0.9 $ a 
5 z, ae . Se s 
E © so "o 3 5 So Ped 
0.1 0.8 
ò o o] a] S ò o a] S 
Kà a? © Y oF Kà ae © Y oP? 
Ss a a a s A a a 
cx e è è X e e e 
Constrained networks generate less Constrained networks have 
preserved latent dynamics different connectivity structure 
e 1.0 f 
144 
10 5 
< 08 8 > 8 124 
= 06 oo FEST 
E E Standard 8 £ 6 5 5 
= 0.4 a=0 E Be eg 
è} Ba=0.05 & Dy 82 
0.2 E a = 0.25 S az 6] 
0 E q = 0.50 2 44 
ee ee es | CT S. 
12 4 6 8 10 9606 >] ESN] 
Neural mode a Y o? vA Soh 9? 
a 4 v “4 4 
eo è & © é ó 


Fig. 5 | Similar behavioural output is necessary but not sufficient for the 
preservation of latent dynamics. a, We trained two sets of RNN models to 
perform the monkey centre-out reaching task. The first set of models were 
standard RNNs, whereas the second were constrained to produce distinct 
latent dynamics from the standard networks. The weight of this extra constraint 
was controlled by parameter a. b, Example movement trajectories for standard 
and constrained networks. c,d, Both standard and constrained networks 
produced similar behavioural output, as shown by the error (c) andthe 
correlation between ‘reach trajectories’ across networks (d). Datashowten 
networks initialized from different random seeds for each type of network 
(n=45 comparisons within each network condition). Error bars, mean +s.d. 

e, Increasing the value of a decreased the preservation of latent dynamics, 
establishing that behavioural similarity is necessary but not sufficient to have 
preserved latent dynamics. Line and shaded area, mean + s.d. f, Networks with 
distinct latent dynamics show differences in underlying connectivity. The 
variance and the dimensionality of the weight changes increased as the latent 
dynamics were constrained to be more distinct from those of standard 
networks by increasing a. Variance values are scaled by 10%. Line and shaded 
area, mean + s.d. 


produce latent dynamics as a solution for behavioural output. In this 
way, preserved latent dynamics across individuals reflect a fundamental 
property of the neural basis of behaviour. 

Here, we introduced two competing possibilities to explain how 
similar behaviour can emerge from different neural circuitry. In the 
first hypothesis, each individual’s brain must learn during develop- 
ment to produce latent dynamics within their idiosyncratic neural 
circuits to enable the desired behaviour, with no guarantee that both 
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individuals arrive at the same ‘solution’. In the second hypothesis, the 
generated latent dynamics for the desired behaviour are constrained 
by the genetically specified organization at the circuit level shared 
by each individual of a given species. Our results support this second 
hypothesis and raise an intriguing possibility: given that in higher ver- 
tebrates the genome does not specify the implemented architecturein 
great detail, for example, to the level of the synapse? *, the genome may 
provide a ‘generative model that is instantiated by each individual’s 
brain. This generative model may constrain low-level details of the 
neural circuitry such that the appropriate neural population latent 
dynamics required for the behavioural repertoire emerge throughout 
development. 
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Methods 


Subjects and behavioural tasks 

Monkeys. We trained four monkeys (monkeys C, MandJ: male, Macaca 
mulatta; monkey T: male, M. fascicularis; aged 6-10 years) to sit in 
a primate chair and make reaching movements using a customized 
planar manipulandum. The movement of a cursor on a computer 
screen was mapped to the motion of the handle of the manipulan- 
dum and the behavioural task was run through custom software in 
Matlab (The Mathworks). Monkeys C, M and J were trained to perform 
atwo-dimensional centre-out reaching task for at least several months 
before the neural recordings, ensuring they had reached expert 
performance. Monkeys C, M and T were trained on a more complex 
random target sequential reaching task. In the centre-out task, the 
monkey moved its hand to the centre of the workspace to begin each 
trial. After a variable waiting period, the monkey was presented with 
one of eight outer targets. The targets were equally spaced ina circle 
and selected randomly with uniform probability. Then, an auditory 
go cue signalled the animals to reach to the target. Monkeys were 
required to reach the target within 1s after the go cue and hold for 
0.5 sto receive aliquid reward, except for monkey J, who was trained 
without the instructed-delay period or the 0.5 s target hold time and 
therefore made larger movements (Extended Data Fig. 1a, right). For 
the centre-out task, there were 12 sessions for monkey C, 6 sessions 
for monkey M and 3 sessions for monkey J. 

In the random target task, the monkeys made four consecutive 
reaches to random targets within a10 x 10 cm? workspace in each trial. 
Each target was presented sequentially ina random location within an 
annulus with 5 cm inner radius and 15 cm outer radius of the previous 
target to enforce minimum and maximum reach lengths. Monkeys 
received a liquid reward during a short break after each successful 
sequence of four random target acquisitions. There was no explicit audi- 
tory go cue and only a brief hold period within the target (100 ms) and 
thena brief delay period (100 ms) before the next target was presented. 
These short constraints helped to enforce that the monkeys made 
separate, directed movements but did not require that the monkey 
necessarily stop between movements. For the random target task, there 
was one ‘reference’ session for monkey C, six sessions for monkey M 
and four sessions for monkey T. As the monkeys performed these tasks, 
we recorded the position of the endpoint at a sampling frequency of 
1kHz using encoders in the joints and digitally logged the timing of 
task events, such as the go cue. Portions of the centre-out reaching 
data have been previously published and analysed in refs. 26,28,46,65. 
Portions of the random target data have been previously published 
and analysed in refs. 31,32. 


Mice. Four 8-16-week-old mice were trained to perform a forelimb 
reaching and pulling task (similar to refs. 38,66) for approximately one 
month, following habituation to head-fixation and the recording setup. 
In each trial, mice had to reach and pull a joystick positioned about 
1.5 cm away from the initial hand position. The joystick appeared, with- 
out any cue, in one of two positions (left or right, less than 1 cm apart). 
Mice could then self-initiate a reach to the joystick and pull it inwards 
to get aliquid reward. Thejoystick was weighted with either a3 ora6g 
load (light or heavy), making up a total of four trial types (two joystick 
positions by two loads). Each trial type was repeated 20 times before 
task parameters were switched to the next trial type without any cue. 
Each session consisted of two repetitions of each set of four trial types 
presented inthe same order, making up 2 x 4 x 20 = 160 trials. Trials with 
incorrect responses (for example, pushing the joystick past a thresh- 
old, 5 mm) or timeout (the lack of pull or push for 10 s) were marked as 
unsuccessful. All joystick operations were programmatically controlled 
using acustom-written open-source Python package: (https://github. 
com/janelia-pypi/mouse_joystick_interface_python). Mice were main- 
tained on a 12/12 h (08:00-20:00) light/dark cycle and recordings were 


made between 09:00 and 15:00. The holding room temperature was 
maintained at 21+ 1°C witha relative humidity of 30-70%. 

There were two sessions for mouse 38, one session for mouse 39, 
two sessions for mouse 40 and one session for mouse 44. Movement 
kinematics were tracked using markerless video-based pose estimation. 
Annotation of behaviour was accomplished using Janelia Automatic 
Animal Behavior Annotator”. Briefly, behaviour was recorded using 
two synchronized high-speed (500 frames s”), high-resolution mono- 
chrome cameras (Point Grey Flea3; 1.3 MP Mono USB3 Vision VITA 1300; 
Point Grey Research) with 6-15 mm (f/1.4) lenses (C-Mount), placed 
perpendicularly in front and to the right of the animal. A custom-made 
near-infrared light-emitting diode light source was mounted on each 
camera. Video was recorded using custom-made software developed 
by the Janelia Research Campus Scientific Computing Department and 
IO Rodeo. This software controlled and synchronized all facets of the 
experiment. For the main analyses, light and heavy trials were pooled 
together because we focused on the reaching phase of the task and 
the location of the joystick does not depend on its weight. Note that 
in Extended Data Fig. 6a we repeated the main analysis to demonstrate 
preserved latent dynamics during the pulling phase, considering all 
four conditions. 


Neural recordings 

Monkeys. All surgical and experimental procedures were approved 
by the Institutional Animal Care and Use Committee of Northwestern 
University under protocol no. ISO0000367. We implanted 96-channel 
Utah electrode arrays in the primary motor cortex (M1) or dorsal pre- 
motor cortex (PMd) using standard surgical procedures. Throughout 
the paper, neural recordings from these two subregions were pooled 
together and denoted as motor cortex. This allowed us to ensure that 
we could evaluate overt and covert dynamics within the same popula- 
tion. Implants were done in the opposite hemisphere of the hand ani- 
mals used in the task. Monkeys M and T received two arrays in M1 and 
PMd simultaneously. Monkey J received a single array in M1. Monkey C 
received two sets of implants: one array in the right M1 while perform- 
ing the task using the left hand and, following removal of this original 
implant, two arrays simultaneously in the left M1 and PMd while using 
the right hand (respectively, monkeys C, and C, in our previous work”). 
Note that for all across-individual analyses, C, and C, are considered 
the same animal. 

Neural activity was recorded during the behaviour using a Cerebus 
system (Blackrock Microsystems). The recordings on each channel 
were band-pass filtered (250 Hz-5 kHz) and then converted to spike 
times on the basis of threshold crossings. The threshold was set to 
5.5x the root-mean-square activity on each channel. We also manually 
spike sorted the recordings from monkeys C, Mand T to identify puta- 
tive single neurons. Monkey J had fewer well-isolated single units than 
the other monkeys, so rather than spike sorting we directly applied the 
multi-unit threshold crossings acquired on each electrode. However, 
it has been shown that the latent dynamics estimated from multi-unit 
and single neuron activity are similar®, an observation that holds true 
for aligning latent dynamics with CCA” (note that we refer to both 
single neurons and multi-units simply as units). We included multiple 
experimental sessions from each monkey: for the centre-out reaching 
task, eight from monkey C,, four from monkey Ca, six from monkey M 
and three from monkey J (example data in Extended Data Fig. 1); for 
the random target task, one ‘reference session’ from monkey C, six 
from monkey M and four from monkey T (example data in Extended 
Data Fig. 8). These experimental sessions were chosen on the basis 
of the high number of units or trials and blind to the behaviour of the 
animal. For the centre-out reaching task, the average number of units 
included for each monkey was: monkey C,, 277 + 14 (mean + s.e.m.; 
range, 210-345); monkey Cp, 85 + 4 (range, 73-92); monkey M, 117 + 4 
(range, 106-130); and monkeyJ, 63 + 9 (range, 54-81). For the random 
target task, the average number of units included was: monkey C,, 280 


(one session only); monkey M, 127 + 9 (range, 101-153); and monkey T, 
49 + 8 (range, 30-66). A more detailed description of the behavioural 
and neural recording methods is presented in ref. 26. 


Mice. All surgical and experimental procedures were approved by 
the Institutional Animal Care and Use Committee of Janelia Research 
Campus. A brief (less than 2 h) surgery was first performed to implant 
a three-dimensional-printed headplate®. Following recovery, the 
water consumption of the mice was restricted to 1.2 ml per day, totrain 
them in the behavioural task. Following training, a small craniotomy 
for acute recording was made at 0.5 mm anterior and 1.7 mm lateral 
relative to bregma in the left hemisphere. A neuropixels probe was 
centred above the craniotomy and lowered with a 10° angle from the 
axis perpendicular to the skull surface at a speed of 0.2 mm min”. The 
tip of the probe was located at 3 mm ventral from the pial surface. 
After a slow and smooth descent, the probe was allowed to sit still at 
the target depth for at least 5 min before initiation of recording to allow 
the electrodes to settle. 

Neural activity was filtered (high-pass at 300 Hz), amplified (200x 
gain), multiplexed and digitized (30 kHz) and recorded using the 
SpikeGLX 3.0 software (https://github.com/billkarsh/SpikeGLX). 
Recorded data were preprocessed using an open-source software 
KiloSort 2.0 (https://github.com/MouseLand/Kilosort) and manu- 
ally curated using Phy (https://github.com/cortex-lab/phy) to identify 
putative single units in each of the primary motor cortex and dorso- 
lateral striatum. A total of six experimental sessions (from four mice; 
Extended Data Fig. 5) with simultaneous motor cortical and striatal 
recordings were included in this work. The average number of motor 
cortical units included for each mouse was: mouse 38, 98 + 4 (range, 
95-102); mouse 39, 64; mouse 40, 75 + 5 (range, 70-80); and mouse 44, 
55. The average number of striatal units included for each mouse was: 
mouse 38, 100 + 13 (range, 87-112); mouse 39, 108; mouse 40, 74 +5 
(range, 69-79); and mouse 44, 110. 


Data analysis 
We used a similar approach for both monkey and mouse data. In all the 
analyses, we only considered the trials in which the animal successfully 
completed the task within the specified time and received a reward. 
We concatenated trials in time for subsequent analyses—that is, no 
trial-averages were taken. For the monkey centre-out reaching task and 
the mouse reaching and pulling task, an equal number of trials to each 
target was randomly selected (eight targets for the monkeys and two 
targets for mice, except in Extended Data Fig. 6a, for which four targets 
were considered). Trial order was randomized to eliminate the possible 
effect of the passage of time. Within each trial, we isolated a window of 
interest that captured most of the movement, starting 50 ms before 
movement onset and ending 400 ms after movement onset. To analyse 
covert behaviour in monkeys, we used a window that spanned the move- 
ment planning period, which started 400 ms before movement onset 
and ended 50 ms after movement onset. Importantly, all of our results 
held when changing the analysis windows within a reasonable range. 
For the monkey random-walk task, each reach could start and end 
anywhere within the workspace. To define movements (conditions) that 
could be matched across animals, we first segmented the workspace 
into 12 circular subsections. Each subsection was then divided into 
six equal sectors and targets in the same angular sector were grouped 
together, creating 72 possible target conditions. We separated the 
sequences of four consecutive reaches and considered each reach as 
a separate movement. To assign each movement toa target condition, 
we first assigned each movement to one of the subsections on the basis 
of the starting position of the given movement, excluding movements 
that started more than 2 cm fromthe centre of the subsection. We then 
recentred the movements so that they started in the centre of each sub- 
section and reached outwards towards their target. The movement was 
then assigned to asector and target condition on the basis of the angle 


of target. To study the preservation of latent dynamics across monkeys 
performing similar behaviour, we needed to match movements (reach 
conditions) across sessions for different monkeys. To maximize the 
number of matched movements, we compared all sessions for Monkey 
Mand Monkey T against a reference session for Monkey C, that had the 
most successful trials. We matched movements in each pair of sessions 
by minimizing the mean squared error (MSE) between pairs of move- 
ments, excluding matches that had MSEs above the threshold of 2% of 
MSEs calculated for all possible pairs of movements. If the matched 
movements had different corresponding target conditions, we used 
the target condition label from the reference session. After this pro- 
cess was completed, we excluded target conditions with less than six 
matched movements, such that paired sessions had up to 29 shared 
target conditions. Because these movements were more ballistic than 
in the centre-out task, we examined a window starting at movement 
onset and ending 350 ms after movement onset. 

All the analyses were implemented in Python using open-source 
packages such as numpy, matplotlib, sci-kit, scipy and pandas” and 
custom code. As we were analysing existing datasets on an individual 
basis, no explicit planning of sample size, group randomization or 
blinding was performed. 


Behavioural correlation. To assess the behavioural stereotypy of a 
given animal, we calculated hand trajectory correlations (Pearson’s 
correlation) of every pair of trials within a session (Extended Data Fig. 1b 
and Extended Data Fig. 5b). The distributions in Fig. 2k inset illustrate 
these correlations pooled across all the monkey centre-out and mouse 
reaching and pulling sessions included in this work. To determine the 
behavioural similarity across pairs of sessions from different monkeys 
or mice (Fig. 2k), we similarly calculated correlations to compare all 
pairs of trials from the two sessions. 


Neural population latent dynamics. To estimate the latent dynamics 
associated with the recorded neural activity in each session for both 
mice and monkeys, we computed a smoothed firing rate as a func- 
tion of time for each unit. We obtained these smoothed firing rates 
by applying a Gaussian kernel (o = 50 ms) to the binned square-root 
transformed spike counts (bin size 30 ms) of each unit. We excluded 
units with alow mean firing rate (less than 1 Hz mean firing rate across 
all bins) but we did not perform any further exclusions, for example, 
based on lack of modulation or behavioural tuning. For each session, 
this produced a neural data matrix X of dimension n by T, where nis the 
number of recorded units and 7 the total number of time points from 
all concatenated trials ona given day; Tis thus given by the number of 
targets per day x number of trials per target x number of time points 
per trial. We performed this concatenation as described above after 
randomly subselecting the same number of trials for all targets for each 
animal (15 trials for monkey centre-out, six for monkey random walk, 
22 for mouse reaching and pulling). For each session, the activity of n 
recorded units was represented as a neural space—an n-dimensional 
sampling of the space defined by the activity of all neurons in that brain 
region. In this space, the joint recorded activity at each time bin is rep- 
resented as a single point, the coordinates of which are determined by 
the firing rate of the corresponding units. Within this space, we esti- 
mated the low-dimensional latent dynamics by applying PCA to X. This 
yielded n PCs, each a linear combination of the smoothed firing rates 
of alln recorded units. These PCs are ranked on the basis of the amount 
of neural variance that they explain. We defined an m-dimensional, 
session-specific manifold by only keeping the leading m PCs, which 
we referred to as neural modes. We chose a manifold dimensionality 
m=10, based on previous studies examining motor cortical recordings 
during upper limb tasks*****. Across all datasets, a ten-dimensional 
manifold explained about 60% of the neural variance for each of the 
monkey motor cortex (Extended Data Fig. 1c), mouse motor cortex 
and mouse striatum (Extended Data Fig. Se). Note, however, that our 
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results held within a reasonable range of dimensionalities, similar to 
refs. 26,33,46 (Extended Data Figs. 2fand 4b). We computed the latent 
dynamics within the manifold by projecting the time-varying smoothed 
firing rates of the recorded neurons onto the m = 10 PCs that span the 
manifold. This produced a data matrix L of dimensions m by T. 


Aligning latent dynamics through CCA. We addressed our hypothesis 
that different animals performing the same behaviour would share 
preserved latent dynamics by aligning the dynamics using CCA”, 
CCA was applied to the latent dynamics of each pair of sessions 
after concatenating the same number of randomly ordered trials to 
each target (condition, in the case of the sequential reaching task). For 
details on using CCA to align latent dynamics, see ref. 26. 

We measured the similarity in latent dynamics across animals by 
computing the across-animal correlations as the canonical correlations 
(CCs) across all pairs of sessions from any two different monkeys or 
mice. To establish the strength of the across-animal correlations, we 
computed an upper bound defined by the within-animal correlations, 
which we calculated as the 99th percentile of the correlations between 
two randomly selected subsets of trials within any given session over 
1,000 samples. The ‘control’ correlations represent a lower bound for 
the CCs. We computed these by shuffling the targets across the two 
sessions and using a randomly selected control window (more details 
in the ‘control analyses’ section below) in each trial, rather than the 
movement or preparatory epochs. 

Note that to summarize each comparison to a single datapoint (for 
example, in Fig. 2k and Extended Data Figs. 2h and 6d), we computed 
the mean of the top four CCs of the latent dynamics”. In Fig. 2k, we 
used this approach to establish a relationship between the strength of 
preservation of the latent dynamics and the consistency of behaviour, 
quantified as the mean trajectory correlation of all possible pairs of 
trials across two animals. Furthermore, when showing pairs of ‘aligned’ 
trajectories across animals, suchas in Fig. 2e and Extended Data Fig. 3, 
the CCA axes were made orthogonal using singular value decomposi- 
tion for visualization purposes. 

Finally, we showed that preserved latent dynamics could be 
uncovered across a broad range of manifold dimensionalities. In 
Extended Data Fig. 2f we repeated the alignment analysis for manifold 
dimensionalities m = 2-19. 


Decoding analysis. To test whether the aligned latent dynamics 
maintain movement-related information, we built standard decod- 
ers to predict hand trajectory during overt behaviour. If the aligned 
latent dynamics across different animals were behaviourally relevant, 
they would allow predicting time-varying hand trajectories even if 
the methods used to identify them (PCA and CCA) are not supervised, 
that is, they do not attempt to optimize decoding performance. We 
compared the predictive accuracy of three different types of decod- 
ers: (1) a within-animal decoder trained and tested (using ten-fold 
cross-validation) on two non-overlapping subsets of trials from each 
session of each animal; (2) an across-animal ‘aligned’ decoder that was 
trained onthe aligned dynamics from one animal and tested on another, 
acomparison we performed on each pair of sessions from two different 
animals; (3) an across-animal ‘unaligned’ decoder that was trained on 
the latent dynamics from one animal and tested on another without 
aligning the dynamics using CCA. We also performed a similar analysis 
to predict the upcoming target during covert movement preparation 
in monkeys (Fig. 4f). 

Hand trajectory decoders were LSTM models with two LSTM lay- 
ers, each with 300 hidden units, followed by a linear output layer. The 
models were implemented with Pytorch” and trained for ten epochs 
with the Adam optimizer, with a learning rate of 0.001. Upcoming tar- 
get classifiers were Gaussian Naive Bayes models” (the GaussianNB 
class in ref. 72). We included three bins of recent latent dynamics his- 
tory, for a total of 90 ms, in the input of both the decoders and the 


classifiers. These extra neural inputs incorporate information about 
intrinsic neural dynamics and account for transmission delays. The R? 
value, defined as the squared correlation coefficient between actual 
and predicted hand trajectories, was used to quantify decoder perfor- 
mance. Moreover, in Extended Data Fig. 4d we verified that our choice of 
across-animal decoder accuracy metric did not influence the observa- 
tion that preserved latent dynamics are informative about behaviour 
by also computing a variance accounted for (VAF) metric, defined as: 


ya (9, -y) 
VAF=1- =! —, 
die (y,-Y) 

where y; represents the actual value of the predicted variable, 9; its 
predicted value and y its mean. For this analysis, we normalized hand 
trajectories by the length of the reaches (determined by the 99th per- 
centile of their hand positions along each axis) because monkeys had 
workspaces of different sizes. 

The hand trajectory was a two-dimensional signal in monkeys anda 
three-dimensional signal in mice. We built separate decoders to predict 
hand trajectories along the x, y (and zfor mice) axes. We then reported 
the average performance across all axes. For target classification, we 
reported the mean accuracy of the classifier (the score() method). 

To test how many dimensions of the aligned latent dynamics were 
needed for accurate across-animal decoding of behaviour, we repeated 
the decoding analysis in the monkey centre-out dataset for manifold 
dimensionalities m =1, 2...,14 (Extended Data Fig. 4b). 

Finally, we performed acontrol analysis to ensure our across-animal 
decoding results were not biased by sharing similar trials for both align- 
ment and decoder training. We split the full dataset of one animal into 
three non-overlapping sets: one to align the latent dynamics, one to 
train the decoder and one to test the performance across animals. 
Extended Data Fig. 4c shows the result of this analysis for the monkey 
centre-out data. Despite having aligned the latent dynamics only using 
half of the data, the impact on decoding performance is negligible. 


Control analyses 

Alignment of latent dynamics with random behavioural windows. 
To establish a ‘behaviourally irrelevant’ window as control, we ran- 
domly selected windows of similar length to our behavioural windows 
(450 ms) along the entire duration of the intertrial and trial periods 
combined. This ensured we had samples of dynamics in the neural 
population with realistic statistics but that they were not directly cou- 
pled to shared behaviour across individuals. We used this window to 
provide alower-bound control for the alignment of neural population 
latent dynamics (‘control in Figs. 2f,g,j, 3d,e and 4b,e and Extended 
Data Figs. 2b-d,g, 3 and 8d). 


Aligning latent dynamics through Procrustes analysis. We used CCA 
to align the latent dynamics in all the analyses. However, to ensure that 
our results hold regardless of the specific method used for alignment, 
we replicated the main result using Procrustes analysis”. Procrustes 
finds the best transformation that minimizes the sum of the squares of 
the differences between the two input datasets. Following a procedure 
identical to the CCA analysis, we aligned the dynamics from two differ- 
ent datasets using Procrustes analysis (the scipy.spatial.procrustes 
class in ref. 73) and then correlated the aligned dynamics to yield a 
metric comparable to that of the CCA (Extended Data Fig. 2g,h). Note 
that the degrees of preservation of latent dynamics obtained with CCA 
and Procrustes analysis are largely similar. 


Neural variance explained by aligned latent dynamics. We meas- 
ured the percentage of neural variance explained by the preserved 
latent dynamics using a method we devised in ref. 33. Briefly, we 
‘reconstructed’ the preserved neural activity by projecting the aligned 
latent dynamics along the CC axes back to the PCA space (the neural 


manifold) and then to the original neural state space. We then measured 
the difference between the total neural variance and the variance of 
these reconstructed signals using an approach similar to that in ref. 78. 
By repeating this procedure iteratively for an increasing number of 
manifold dimensions m, we measured the neural variance explained 
by each dimension of the aligned latent dynamics. Using this approach, 
we found that preserved latent dynamics explain a significant fraction 
of the neural population variance (Extended Data Fig. 2e). 


Surrogate datasets with TME. We established alower-bound control by 
aligning the latent dynamics from randomly selected windows sampled 
across different task conditions and behavioural epochs (see above). In 
addition to this control, we also used TME to generate surrogate neural 
data as another lower-bound control”. TME produces surrogate data 
that preserve the second-order statistics of the actual neural data (that 
is, covariance across time, across neurons or across experimental condi- 
tions) but are otherwise random (Extended Data Fig. 2a). Aligning these 
surrogate data through the same procedure as the original data shows 
significantly lower correlations for monkey centre-out task, random-walk 
task and mouse reaching and pulling task (Extended Data Fig. 2b-d). 


Aligning topological structure in neural population activity. To test 
whether the topological structure in the produced movements is suffi- 
cient to produce preserved latent dynamics, we quantified the degree of 
similarity in latent dynamics across individuals that could be uncovered 
whenaligning the static, topological features of the neural population 
activity, rather than the dynamics of the movements, using a technique 
developed in ref. 26. To align the topological structure of neural popula- 
tion activity, we time-averaged the activity for each neuron during the 
execution epoch of each trial in the monkey centre-out reaching task. 
We then analysed the time-averaged data with the previous methodol- 
ogy by performing PCA to find a neural manifold and using CCA to align 
each pair of sessions (Extended Data Fig. 7a). This procedure led to 
well-aligned ‘topological representations’ (example in Extended Data 
Fig. 7b). To directly test whether aligning the topological structure of 
neural population activity is sufficient to uncover preserved latent 
dynamics, we projected the latent dynamics on the CC axes found 
through this (static) topological alignment and calculated the pairwise 
correlations of the resultant projected latent dynamics. These correla- 
tions were remarkably lower than those obtained through alignment of 
the time-varying latent dynamics (Extended Data Fig. 7c,d). 


Control analyses on the numbers of conditions and neurons. To 
establish that the preservation of latent dynamics holds across different 
degrees of task complexity, we calculated the correlations for increas- 
ing numbers of subsampled target conditions for each pair of sessions 
inthe monkey random target task (Fig. 3f and Extended Data Fig. 8b). 
We randomly subsampled different combinations of target conditions 
and calculated the degree of preservation of the latent dynamics for up 
to 10,000 combinations for each number of conditions. 

To establish that preserved latent dynamics can be uncovered regard- 
less of the specific measured neurons, we also calculated the correlations 
for varying numbers of neurons in the random target task (Fig. 3g and 
Extended Data Fig. 8c). For each pair of sessions, we either randomly kept 
neurons (Fig. 3d) or randomly dropped neurons (Extended Data Fig. 8c) in 
increments of ten until we ran out of measured neurons for either session 
and repeated this process 50 times, calculating the degree of preserva- 
tion at each step. For both analyses, we calculated the mean correlations 
for the top four CCs across all subsamples for each pair of sessions. 


Comparison of different but related tasks. The central hypothesis of 
this study is that preserved latent dynamics are the basis for the gen- 
eration of similar behaviour across individuals from the same species. 
Here, we sought to further support this hypothesis by showing that 
the latent dynamics produced by two individuals engaged in the same 


task are more similar than the latent dynamics produced by the same 
individual performing two different but related tasks. To this end, we 
compared our results to our previous study onthe relationship of neural 
population activity underlying different but related wrist manipulation 
or reach-to-grasp tasks in monkeys” (Extended Data Fig. 9). 


Recurrent neural network models 

Model architecture. To show that the preservation of latent dynamics 
across animals engaged in the same task is not a trivial consequence 
of similar behaviour, we trained RNNs to perform the same centre-out 
reaching task as the monkeys. These models were implemented using 
Pytorch”. Similar to previous studies simulating motor cortical dynam- 
ics during reaching?”” *!, we implemented the dynamical system 
x =F (x, s) to describe the RNN dynamics: 


N I 
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where x; is the hidden state of the ith unit and r, is the corresponding 
firing rate following tanh activation of x;. All networks had N = 300 units 
and /=3 inputs, a time constant t= 0.05 s and an integration time step 
dt=0.01s. The noise n was randomly sampled from the Gaussian distri- 
bution .\(0,0,2)for each time step. Each unit had an offset bias, b;, which 
was initially set to zero. The initial states X,- were sampled from the uni- 
form random distribution 2/(—0.2,0.2) . All networks were fully recur- 
rently connected, with the recurrent weights/ initially sampled fromthe 
Gaussian distribution n(o, +) , where g= 1.2. The time-dependent 


inputs s fed into the network had input weights B initially sampled from 
the uniform distribution ¢/(-0.1,0.1) . These inputs consisted of a 
one-dimensional fixation signal which started at 2 and went to O at the 
go cue and a target signal that remained at O until the visual cue was 
presented. The two-dimensional target signal (2 cos 0°"®", 2 sin 0"8"') 
specified the reaching direction 0"*" of the target. 

The networks were trained to produce two-dimensional outputs 
p corresponding to x and y positions of the experimentally recorded 
reach trajectories, which were read-out via the linear mapping: 


N 
PAt)= È Wer (0) 
kal 


where the output weights Wwere sampled from the uniform distribu- 
tion U(-0.1,0.1). 


Model training. Networks were trained to generate positions of reach 
trajectories using the Adam optimizer® witha learning rate / = 0.0005, 
first moment estimates decay rate £, = 0.9, second moment estimates 
decay rate £, = 0.999 and € = 1 x 108. The loss function was defined as 
the MSE between the two-dimensional output and the target positions 
over each time stept, with the total number of time steps T = 400. The 
first 50 time steps were not included to allow network dynamics to relax: 


1 B T 
L= — 2 2 2 (ptarset (b,t) — poutput (b, t))?. 
2B(T-50) 6H) £50 d=1.2 i d 


To examine whether two networks could have different latent dynam- 
ics while producing the same motor output, we devised a network 
with more constraints to perform the behavioural task with distinct 
latent dynamics (Fig. 5a). We added a loss term that penalized the CC 
between the latent dynamics of the ‘constrained’ network being trained 
and those of another previously trained ‘standard’ network during 
movement execution: 
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where c,is the ith CC. To examine different degrees of preserved latent 
dynamics, we trained the networks at varying values of a = 0, 0.05, 
0.25 or 0.50. 

Networks were trained until the average loss of the last ten training 
trials fell below a threshold of 0.2 cm’, for at least 50 and up to 500 
training trials, with a batch size B = 64. Each batch had equal numbers 
of trials for each reach direction. We clipped the gradient norm at 0.2 
before the optimization step. Both standard and constrained training 
were performed on ten different networks initialized from different 
random seeds. The same set of random seeds was used for constrained 
networks at different values of a. 


Connectivity analyses. By increasing the value of a, we were able to 
decrease the preservation of the latent dynamics while keeping behav- 
ioural performance constant. To examine how this changed the underly- 
ing connectivity, we calculated the variance and dimensionality of the 
weight changes in the recurrent weights/ following training (Fig. 5f,g). 


Statistics and reproducibility 

We compared the performance of various within-animal and across- 
animal movement decoders and classifiers using two-sided Wilcoxon’s 
rank sum tests. We replicated the core findings across two species (mice 
and monkeys), four behaviours (a centre-out reaching task, a sequen- 
tial reaching task and a reach, grasp and pull task, along with during 
covert movement planning) and two brain regions (motor cortex and 
dorsolateral striatum). Experiments on each species were performed 
independently in two different laboratories and by different scien- 
tists. The mice experiments were done ina single cohort, whereas the 
monkey data were collected in two sets of experiments (one for the 
centre-out task, another for the random reaching task), each spanning 
2 years. Overall, our neural recordings and behavioural data are in good 
agreement with related published studies. All attempts at replication 
were successful. 


Reporting summary 


Further information on research design is available in the Nature 
Portfolio Reporting Summary linked to this article. 


Data availability 


Most of the monkey datasets are publicly available on Dryad (https:// 
datadryad.org/stash/dataset/doi:10.5061/dryad.xd2547dkt) and 
CRCNS (https://doi.org/10.6080/KOFT8J72). The remaining monkey 
datasets and the mouse datasets will be made available on reasonable 
request. 


Code availability 


All analyses were implemented using custom Python code and using 
open-source software. All the result panels are reproducible by running 
Jupyter notebooks. Code to reproduce all the results is openly avail- 
able in https://github.com/BeNeuroLab/2022-preserved-dynamics. 
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Behavioural similarity across monkeys 
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Extended Data Fig. 1| Behavioural and neural data from monkeys 
performing acentre-out reaching task. A) Hand trajectories ofa single 
session of every monkey included in this work. Each line represents one 

trial, colour coded by target (inset). Scale bar, 1 cm. B) Distribution of the 
behavioural correlations for all pairs of trials from every included session 

for each of the three monkeys. Grey circle, mean. Blue triangle, mean for the 
representative session shown in (A). C) Cumulative neural variance explained 
as a function of the number of neural modes included. Data include 12 sessions 
for monkey C, 6 sessions for monkey M and 3 sessions for monkeyJ. Each line, 
one session. Dashed blue line, the same session shown in (A). 
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Preserved latent dynamics exceed a Tensor Maximum Entropy (TME) control that replicates 
the statistics of the recorded neural data 
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Preserved latent dynamics across individuals for a 
range of neural manifold dimensionalities 
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Preserved latent dynamics can be uncovered with other linear alignment methods 
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Extended Data Fig. 2 | Additional verification of preserved latent 
dynamics. A) We used Tensor Maximum Entropy (TME) to generate surrogate 
data that preserved all the second-order statistics of the actual neural data. 
Element-wise conservation of the covariance across time, targets and neurons 
between recorded and surrogate data for one example monkey centre-out 
dataset. B) Preserved latent dynamics across all pairs of sessions from three 
different monkeys performing the centre-out reaching task. Data include all 21 
sessions and 126 comparisons shown in Fig. 2g. Note that the correlation 
between aligned latent dynamics across monkeys (red) was quite close to the 
within-monkey correlations (grey) and largely exceeded lower-bound controls 
from both randomly sampled behavioural epochs (orange) and TME (mustard). 
Line and shaded area, mean + s.d. C) Similar to (B) but for mouse motor cortex 
during the reaching and pulling task. Data include all six sessions and 13 
comparisons across four different mice shown in Fig. 2j. Line and shaded area, 
mean +s.d. D) Similar to (B) but for the monkey motor cortex during the 
sequential reaching task. Data include all 11 sessions and ten comparisons 
across three different monkeys shown in Fig. 3e. Line and shaded area, mean + 
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s.d. E) Preserved latent dynamics capture dominant features of neural 
population activity. Comparison of the total neural variance explained by 

the top four PCs (neural modes) and the top four Canonical Correlations. 

Single dots, individual comparisons for all monkey centre-out datasets. 

F) Preservation of latent dynamics does not depend on the number of manifold 
dimensions. Average canonical correlation as function of the dimensionality of 
the neural manifold (number of principal components) for four of the datasets 
considered (legend). Data pooled across all comparisons separately for each 
dataset. Line and shaded area, mean + s.d. across the top four canonical 
correlations. MCx, motor cortex. Line and shaded area, mean +s.d.G) Preserved 
latent dynamics can be uncovered witha different linear alignment method, 
Procrustes analysis (details in Methods). Data presented as in Fig. 2g but using 
Procrustes analysis rather than CCA to align the latent dynamics. Line and 
shaded area, mean + s.d. H) Pairwise comparison between CCA and Procrustes 
analysis. Single dots, individual comparisons. Data includes 21 sessions and 126 
comparisons across three monkeys. 


Additional examples: Preserved motor cortical dynamics during reaching 


Example: Aligned motor cortical 
dynamics during grasping and 
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Extended Data Fig. 3 | Additional examples of canonical correlations and 
preserved latent dynamics between pairs of sessions from two different 
monkeys and mice. A) Example canonical correlations between motor cortical 
latent dynamics for monkeys during movement execution for the centre-out 
reaching task. Shownare five of the individual comparisons included inthe 
pooled results in Fig. 2g. B) Example canonical correlations between motor 
cortical latent dynamics for mice during movement execution. Shown are five 


of the individual comparisons included in the pooled results in Fig. 2j.C) Same 
as (B) for mouse striatum. D) Same as (A) but during movement planning. 

E) Trajectories described by the motor cortical latent dynamics before (left) and 
after (right) alignment for a representative pair of sessions from two different 
mice. Individual lines, mean latent trajectory to each target (colour coded asin 
the inset). F) Similar to (E), for striatal latent dynamics in mice. G) Similar to (E), for 
motor cortical latent dynamics during covert movement preparation in monkeys. 
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Example hand position predictions from monkey motor cortex 
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Extended Data Fig. 4 | Preserved latent dynamics across monkeys are 
behaviourally relevant. A) Example kinematic predictions across animals. 
Actual and decoded hand trajectory during one reach to each of the eight 
targets along the X (horizontal, top) and Y (vertical, bottom) directions; arrows 
ontop of eachcolumn indicate target direction. The figure shows predictions 
using three models along with their respective accuracy (legend). Data from 
the same representative session shown in Fig. 2b. Note that a decoder trained 
ona different monkey after alignment of their latent dynamics (‘Across 
aligned’) performed virtually as well as a decoder trained and tested onthe 
same session from the monkey being tested (‘Within’). B) Decoding of hand 
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position based onthe latent dynamics for manifolds with increasing 
dimensionality. Line and shaded area, mean + s.d. C) Across-animal decoding 
accuracy is not significantly impacted by a ‘full’ cross-validation procedure in 
which we used three sets of trials for training and testing (one for alignment, 
one for decoder training and one for testing). Error bars, mean + s.d. D) Aligned 
latent dynamics also allow for accurate across-animal decoding of movement 
kinematics according to a variance explained metric (variance accounted for; 
Methods). Data points, individual comparisons between two sessions from 
different monkeys (n = 126). Error bars, mean + s.d. 
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Extended Data Fig. 5 | Behaviour and neural activity for mice performing for two different mice (top) and corresponding hand trajectories (bottom). 
a grasping and pulling task. A) Hand trajectories of every session included in Each column, one reach to each of the two directions. Arrows on top of each 
this work (notice that two mice have two sessions each). Each line represents column, reach direction, colour coded as in (A). D) Similar to (C) but showing 
onetrial, colour coded by target (left, yellow; green, right; see inset). Circle, the striatal activity during the same example movements. E) Cumulative neural 
initial hand position. Scale bar, 1 mm. B) Distribution of the behavioural variance explained as a function of the number of neural modes included, 
correlations for all pairs of trials from every session. Grey circle, mean. colour coded by brain region. 


C) Example normalized motor cortical firing rates aligned to movement onset 
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Preserved motor cortical latent dynamics 
across four conditions in mouse motor cortex 


Decoding hand position across mice based on 
their preserved motor cortical dynamics 
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The accuracy of across-animal movement decoding 


depends on the similarity of their behaviour 


g i 
E 
> e 
8 e ” 
g 0-6 oe ° ° 
g cone 
$ 04 2%, ° 
5 ŝ ® Monkey motor cortex 
D e © Mouse motor cortex 
a 0.2 6° e Mouse striatum 
i 7 T F T ¥ T F T bi 1 
0.70 0.75 0.80 0.85 0.90 0.95 


Behavioural correlation 


Extended Data Fig. 6| Further exploration of preserved latent dynamicsin 
motor cortex and striatum. A) Preservation of motor cortical latent dynamics 
during the four conditions (two reach directions x two loads) of the peri-grasp 
interval of the mouse grasping and pulling task (red) extends the results during 
the two reaching conditions (blue) from Fig. 2j. The lower-bound control 
(orange) and the within-animal upper bound (grey) from Fig. 2j are shown here 
for reference. Data include six sessions across four different mice (n=13 
comparisons). Line and shaded area, mean + s.d. B) Decoders trained onthe 
aligned latent dynamics from one mouse allowed predicting hand position 
froma different mouse with reasonable accuracy (blue), considerably better 
than decoders based on the unaligned latent dynamics (green). Lines, 
individual comparisons between two sessions from different mice (n = 6 
sessions and n= 13 pairs), colour coded based onthe correlation in hand 
kinematics across them (legend). Error bars, mean + s.d. C-E) Motor cortical 


and striatal latent dynamics are similarly preserved during the mouse pulling 
and grasping task and allow for comparable decoding of movement kinematics. 
Theacross-animal motor cortical (yellow) and striatal (orange) correlations 
averaged across all comparisons (reproduced from Fig. 2j and Fig. 4b) are 
largely overlapping (C). Line and shaded area, mean + s.d. Pairwise comparison 
of the degree of preservation of the latent dynamics (D) and across-animal 
decoder accuracy (E) for motor cortex and striatum. Single dots, pairs of 
sessions from two different mice (n = 6 sessions and n = 13 pairs). F) Decoding 
performance based on the aligned latent dynamics for each pair of sessions 

for the three stereotypical movement execution datasets (legend) as function 
of their mean behavioural correlation. Similar to Fig. 2k. Dots, individual 
comparisons across pairs of sessions that includes 21 sessions and 126 
comparisons across three monkeys and six sessions and 13 comparisons across 
four mice. 


Aligning the topological structure of neural population activity 
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Alignment based on the topological organisation of neural population activity 
is not sufficient to uncover preserved latent dynamics 
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Extended Data Fig. 7 | Aligning the topological organization of the neural 
datais not sufficient to uncover preserved latent dynamics across 
individuals. A) To align the latent dynamics across sessions based on the 
topological organization of the neural data, we performed CCA onthe 
population activity averaged across all time points on each trial. This process 
created a distribution of points with clear topological structure corresponding 
to the different reaching directions. B) This method provided excellent 
alignment of the topology ona single trial basis. Before alignment (left), the 
projection of each trial onto the neural manifold from each monkey (shown 
here overlaid ontop of each other) looked different. After alignment, similar 
target-specific structure is present across both animals. Individual markers, 
areachtoone of the eight targets (colour code in inset) for Monkey C, (closed 
circles) and Monkey M (open squares). C) Pairwise comparisons of the 
correlations after projecting the latent dynamics onto the CC axes found by 
aligning the topology (vertical axis) and onto the CC axes found by aligning the 
latent dynamics (horizontal axis). Data shown for the top six neural modes (see 
legend for colour code). Each dot represents one comparison; data include 21 
sessions and 126 pairs of comparisons from three different monkeys during the 
centre-out reaching task. All dots lie well below the diagonal (dashed grey), 
indicating that aligning the latent dynamics based on the topology does not 
reach the correlation values obtained by aligning the latent dynamics. 

D) Correlation values were significantly lower when the alignment was based 
on preserving the topological structure of the neural population activity rather 
thanits latent dynamics, illustrating the importance of the precise temporal 
dynamics for uncovering preserved latent dynamics. Same data as (C). 
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Extended Data Fig. 8 | Behaviour, neural data and preservation of latent 
dynamics for monkeys performing a sequential reaching task. A) Example 
normalized neural firing rates aligned to movement onset for two different 
monkeys (top) and corresponding hand trajectories (bottom). Each column, 
one reach to each of the 28 conditions, which we defined based on the starting 
location, duration and direction of the reaches. B) Canonical correlations as 
function of the number of conditions considered (legend). Note that canonical 
correlations stabilize after ~9 conditions are considered. Data fromthe same 


z 0.8-4 117 0.87 i 
] 97 
16 = 06-4 7 Lo, 
10 0.4- of 0.44 
A ] 37 
0.24 az 0.24 
04 04 
T T x T F T X T r 1 
12 10 
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example dataset as in (A). C) Canonical correlations as function of the number 
of neurons (legend). Note that the traces line up well after ~57 neurons are 
considered. Data from the same example dataset as in (A). D) Permuting the 
condition labels (blue) considerably decreases the alignment of latent 
dynamics across animals, bringing the canonical correlations closer to the 
lower-bound control. Data from all 11 sessions and ten comparisons across 
three different monkeys. Line and shaded area, mean +s.d. 


Monkeys performed three tasks (Movement, Elastic and Isometric) that 
required generating wrist torques to acquire similar targets 
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Extended Data Fig. 9 | The latent dynamics of two different monkeys 
performing the same task are more preserved than those of the same 
monkey performing two related tasks. A) We compared the latent dynamics 
produced by monkeys performing three related visually-guided one-dimensional 
wrist manipulation tasks using the same manipulandum. We also compared 

two reaching and grasping tasks (details in ref. 33). B) The organization of 

the targets and the produced cursor trajectories were the same across all three 
wrist manipulation tasks (note the overlapping cursor trajectories). C-E) Here, 
we show behavioural and neural data from the two wrist manipulation tasks 
that had the most similar kinematics ((C); lines colour coded by targets according 
to (B)): the one-dimensional movement task and the one-dimensional spring- 
loaded (‘elastic’) movement task. These two tasks engaged the same muscles 
(D) and the activity of single units was also relatively similar across them 

((E) shows three example units). F) The covariance patterns elicited by these 
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tasks as well as the rest of wrist and reaching tasks was also quite similar across 
behaviours (this figure shows the number of neural modes out of twelve that 
are well-aligned across tasks according toa principal angle analysis). G) Despite 
the similarities in task characteristics (B), motor output (C), engaged muscles 
(D) and even single unit tuning (E), the latent dynamics generated by the same 
monkey as it performed different wrist tasks (mustard and light grey) or 
different reaching task (dark grey) within the same day were less similar than 
those of two different monkeys engaged in the same task (red). The across 
monkey results (red) illustrate preservation of the latent dynamics as monkeys 
performed the same centre-out reaching task; the results are the sameasin 
Fig. 2e only recalculated with 12 components to match the across-task analysis 
in Gallego et al.” Lines and shaded area, mean +s.d. across the comparisons 
indicated in the legend. (A-F) are adapted from ref. 33; the across-task 
alignment data in (G) partly reproduces Supplementary Fig. 7 therein. 
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Extended Data Fig. 10 | Additional data on recurrent neural network 
models. A) Behavioural output and network dynamics for each reach direction 
for example Standard (left) and Constrained (right) networks that were 
penalized for canonical correlation similarity to the latent dynamics of 
Standard networks. Results presented as in Fig. 2b. Data for constrained 
networks represents results for a= 0.50. B) Behavioural output produced by 
the two types of RNNs. The Standard (top) and Constrained (bottom) networks 
produce similar output. C) Behavioural correlation within and across pairs of 
Standard and Constrained networks. Data show 10 networks initialized from 
different random seeds for each type of network, 45 comparisons within 
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network types, 100 comparisons across network types. Markers, individual 
comparisons. Error bars, mean +s.d. D) Similarity in latent dynamics between 
pairs of Standard networks (blue) and Constrained networks (purple) and 
across network types (grey). Note that Constrained networks have correlations 
that are much more dissimilar than Standard networks. Likewise, the correlation 
between the latent dynamics of pairs Standard and Constrained networks are 
much lower than those between pairs of Standard networks. Thus, behavioural 
similarity is not sufficient to obtain highly correlated latent dynamics. Line and 
shaded area, mean +s.d. 
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The majority of the monkey datasets are publicly available on Dryad (https://datadryad.org/stash/dataset/doi:10.5061/dryad.xd2547dkt) and CRCNS (http:// 
dx.doi.org/10.6080/KOFT8J72). The remaining monkey datasets and the mouse datasets will be made available upon reasonable request (these data have yet not 
been appropriately curated for public distribution). 
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We analyzed data from four monkeys (including two different behavioral tasks) and four mice to show reproducibility. Overall, our results 
include as many as 32 sessions across the two monkey datasets, and six sessions for the mice dataset. This sample size is larger than most 
studies in the field —in non-human primate research, a few sessions across n = 2 individuals is considered to be sufficient. In both species, we 
simultaneously recorded from enough neurons (or neural units in the case of one monkey), to reliably estimate the neural population latent 
dynamics that our results are based upon (see, e.g., Churchland et al Neuron 2010, Gallego et al Nature Comms, 2018, Sadtler et al Nature 
2014, Trautman et al Neuron 2019). 


No data were excluded from analysis. 


We replicated the core findings across two species (mice and non-human primates), four behaviors (a centre-out reaching task, a sequential 
reaching task, and a reach, grasp and pull task, along with during covert movement planning), and two brain regions (motor cortex and 
dorsolateral striatum). Experiments on each species were performed independently in two different laboratories, and by different scientists. 
The mice experiments were done in a single cohort, whereas the monkey data was collected in two sets of experiments (one for the centre- 
out task, another for the random reaching task), each spanning two years. Overall, our neural recordings and behavioural data are in good 
agreement with related published studies. All attempts at replication were successful. 


Covariates were not experimentally controlled so randomization was not necessary for our study (because all the animals were assigned to 
the same group). 


Blinding was not necessary for our study as we were not assessing the effect of any manipulation or comparing across conditions. 
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Materials & experimental systems Methods 

n/a | Involved in the study n/a | Involved in the study 
Antibodies [| ChIP-seq 
Eukaryotic cell lines [| Flow cytometry 
Palaeontology and archaeology O MRI-based neuroimaging 


Animals and other organisms 
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Animals and other research organisms 


Policy information about studies involving animals; ARRIVE guidelines recommended for reporting animal research, and Sex and Gender in 
Research 


Laboratory animals We recorded from four non-human primates (three macaca mulatta, one macaca fascicularis; all male, aged 6-10 years) and four 
mice (C57BL/6; aged 8-16 weeks). Mice were maintained on a 12-12 h (8am-8pm) light/dark cycle and recordings were made 
between 9am and 3pm. The holding room temperature was maintained at 21+/-1°C with a relative humidity of 30% to 70%. 


Wild animals No wild animals were used in this study. 

Reporting on sex Sex was not considered in this study. 

Field-collected samples This study did not involve samples collected from the field. 

Ethics oversight All surgical and experimental procedures for the monkey dataset were approved by the Institutional Animal Care and Use Committee 


of Northwestern University under protocol \#IS00000367. All surgical and experimental procedures for the mouse dataset were 
approved by the Institutional Animal Care and Use Committee of Janelia Research Campus. 


Note that full information on the approval of the study protocol must also be provided in the manuscript. 


